Novel potential energy surface-based quantum dynamics of ion–molecule reaction
School of Physics and Electronics, Shandong Normal University, Jinan 250358, China
† Corresponding author. E-mail:
qtmeng@sdnu.edu.cn
1. IntroductionThe ion–molecule reaction is of fundamental interest to the interstellar media, plasmas, planetary ionospheres and combustion processes.[1,2] As an ion–diatomic molecule reaction,
system has been playing an increasingly significant role in both the planetary ionospheres and the earth ionosphere. It can also be served as a typical three-atom reaction system in the field of ion–molecule reaction dynamics. So the
reaction and its isotopic variant reaction
have been investigated extensively both experimentally and theoretically.[3–10]
Experimentally, Fehsenfeld et al.[3] have measured the rate constant (300 K) for the reactions in which the positive ion such as
,
, Ar+, etc. extracts a hydrogen atom from H2 in a pulsed-discharge, flowing-afterglow reaction tube. Compared with the ion–molecule reactions performed by Smith et al.[4] with the similar species at the same temperature using the selected ion flow tube techniques, the rate constant obtained is a little bit big. While for the isotopic variant reaction
, Burley et al.[5] used the ion beam mass spectrometry to examine the kinetic energy, the total reaction cross sections and the intra-molecular isotope branching ratio, and the results are in excellent agreement with those predicted by the Langevin–Gioumousis–Stevenson model.[6]
Theoretically, Martínez et al.[7] constructed a global potential energy surface (PES) of the H2O+(
) system, and using the quasi-classical trajectory (QCT)[8,9] method they obtained the integral cross sections (ICSs) and rate constants at temperatures in a range of 200 K–1400 K for the
reaction, which are in excellent agreement with the experimental results cited from Burley et al.ʼs paper.[5] In 2006, Martínez et al.[10] calculated the ICS of the
reaction by the time-dependent real wave packet method within the helicity decoupling approximation,[11] which reinforces the adequacy of the ab initio analytical PES used. Bulut et al.[12] have obtained the ICS and the rate constant of state-to-state quantum dynamics of reaction
by the time-dependent wave packet-Coriolis coupling (TDWP-CC) calculations.
In 2015, Song et al.[13] constructed a novel electronic ground-state PES of H2O+(
). According to the scanning of this PES, the minimum energy path (MEP) of
is plotted in Fig. 1. According to this PES,[13] we have performed the dynamic calculation of relevant isotopic variant reaction by using the Chebyshev wave packet propagation method.[14–22] The related calculation results are presented in this paper. The rest of this paper is organized as follows. In Section 2, we present the theoretical method for obtaining the quantum dynamic information. In Section 3, we provide and discuss the computing results. Finally, we draw some conclusions from the present study in Section 4.
2. Theoretical methodIn our calculations, we adopt the reactant Jacobi coordinates (R, r, γ), where R is the distance between the
atom and the center of mass of the D2 molecule, r is the D–D distance, and γ is the Jacobi angle between R and r. We express the Hamiltonian of this system as
where the
μR and
μr are the reduced mass relating to the two radial Jacobi coordinates
R and
r respectively,
represents the potential energy of the ground electronic state,
J denotes the total angular momentum, and the squared orbital angular momentum operator is expressed as
Here,
and
are the projections of
and
onto the body-fixed
z axis, respectively,
and
are the corresponding raising (lowering) operators, which are responsible for the Coriolis coupling.
[23]In the present study, we use the modified Chebyshev recursion[15] to propagate the initial wave packet
with
and
being the initial wave packets, which can be chosen as
[18]where
,
, and
δ are the central momentum, the central position, and the width of the initial wave packet, respectively;
N the normalization constant;
vi,
ji, and
li are the initial vibrational, rotational, and orbital angular momentum quantum numbers, respectively;
is the vibrational eigenfunction of the reactant D
2. When dealing with the scattering problems by the time dependent quantum wave packet method, we always choose a damping function to enforce the outgoing boundary condition
[14,15] at the grid edges. In our calculations, the absorption potential is defined as
where
Rd is the starting position of the absorption potential and
dR denotes the intensity of the absorption potential. The dynamic information is obtained just by the flux operator method, i.e., a surface is set to separate the channel of the product from the channel of reactant artificially, and the Fourier transform is used to analyze all the wave packets passing through the surface and then entering into the asymptotic region of the product channel to obtain the product distribution.
After the wave packet evolving with time is obtained, the related dynamic information can be acquired. To avoid calculating the S-matrix elements,[24,25] the total reaction probability is calculated using the following flux-based formula for the Chebyshev propagation:[26,27]
where
defines the dividing surface, which is the position of the flux calculated and
θ is the Chebyshev angle. In the above equation,
is the amplitude of the initial wave packet, and its value is given by
[28]where
is the spherical Hankel function of the second kind,
is given by
with helicity quantum number
being the projection of the total angular momentum
J on the body-fixed
z axis in the Coriolis coupling calculation. In the reaction probability calculations, the Coriolis coupling effect is fully considered in a range of
. For
, to reduce the computing time, we cut the maximum helicity quantum number to
. This value of
provides convergent results of the reaction probabilities for the considered collision energy range.
The ICS can be given by summing the initial state-specified reaction probabilities for all Jʼs partial waves,
in which
ki is equal to
with
being the collision energy.
3. Results and discussionIn this work, the extensive convergence tests have been performed for optimizing the numerical parameters. The obtained parameters are summarized in Table 1. We employ 300 and 200 equidistant grids for the R coordinate in a range of 0.01 a.u.–22.0 a.u. (The unit a.u. is short for atomic unit) and for the r coordinate in a range of 0.5 a.u.–14.0 a.u., respectively. For the Jacobi angle γ, 50 Gauss–Legendre quadrature points are taken between
and π. In order to minimize the spectral range, both the potential and the kinetic energy[29] are truncated at 0.9 hartree (
). The center of the initial packet is set to be
a.u. with an average translational energy of 0.14 eV, and a width of 0.10 a.u. The starting position of the absorption potential is set to be
and
with absorption intensites
and
, respectively. The position of flux calculated is placed at
. As is well known, the propagation steps are very important for a convergent dynamic calculation. The convergence tests show that in this study the number of propagation steps in the range of collision energy up to 1.0 eV is 40000.
Table 1.
Table 1.
 | Table 1.
Parameters of the quantum wave packet calculations (in atomic units, unless otherwise stated).
. |
3.1. Reaction probabilitiesThe reaction probabilities of
reaction for five partial waves (J = 0, 20, 40, 60, 100) are displayed in Fig. 2. From this figure, we can observe the phenomena as follows. Firstly, the reaction probability of J=0 has no threshold energy, which is in agreement with the exothermic reaction. Secondly, because the potential well can induce the complex to form, all of the probabilities show oscillatory features in the whole range of collision energy, and the intensity of the oscillation is very strong in the low energy region, while it becomes less pronounced in the high energy region. Thirdly, comparing the titled reaction with the
reaction,[10] we find that because the mass of D2 is greater than that of H2, the reaction probabilities of the titled reaction are obviously smaller than those of the
reaction for all considered partial waves. Finally, we can find that, as a result of the centrifugal barrier increasing, the reaction threshold increases to the higher energy with the increase of J.
The vibrational excitation plays a significant role in a reaction, no matter whether the reaction is exothermic or endothermic. Figure 3 shows the reaction probabilities of the different vibrational states of the reactant D2 for the titled reaction. As can be seen, the threshold energy of the reaction probability is equal to zero. Similarly because of the existence of the potential well in the reaction path, the oscillatory feature is displayed in the whole collision energy range for the reaction probabilities of different vibrational states. With the increase of the collision energy, the reaction curves of all the vibrational quantum numbers (v = 0, 1, 2) gradually decline. The vibrational states with v = 1, 2 for the reagent have a slight effect on the reaction, which is a little lower than ground vibrational state (v = 0) at low collision energies (
), which means that the D2 vibration does not enhance the reaction probability at low collision energies. The result is consistent with a previous study by Bulut et al.[12] However, we find that at the high collision energies (
), the reaction probabilities of the vibrational states with v = 1, 2 obviously lie above that with v = 0 as shown in Fig. 3, which suggests that the raising of vibrational quantum number can accelerate the reaction.
For the sake of getting the effects of the rotational excitation on the reaction, the reaction probabilities of the
(
) reaction for j = 0, 1, 2 are plotted in Fig. 4. Obviously the reaction probabilities gradually decline with the collision energy increasing. From the three curves of j = 0, 1 and 2, it can be seen that the reaction probabilities of j = 0 lie above the others in the whole energy range due to the increase of the centrifugal barrier. Besides, we can find that the decrease of amplitude takes on a nonlinear relationship with the increase of j. This phenomenon indicates that the rotational excitation of the reagent D2 can impede the reaction, which is in excellent agreement with the conclusion obtained by Gao et al.[30] So we can come to the conclusion that the D2 rotation will prevent the
atom from coming close to the D2 molecule for the exothermic and barrierless reaction.
3.2. Integral cross sectionIn order to test the convergence of this calculation for the titled reaction, the J-dependent partial wave contributions (weighted over a 2J +1 factor) to the ICS at the different collision energies are depicted in Fig. 5. We can see that for all curves of Ec = 0.25 eV, 0.5 eV, 0.75 eV, and 1.0 eV, the partial wave contribution increases with J at first. When the partial wave contribution rises to a highest value, it will decrease sharply with J increasing because of the increasing of the centrifugal potential. The oscillatory features of the curves are in agreement with the reaction probabilities of all partial waves. From the curves, it is clear that higher collision energy implies more partial waves contributing to the ICS. It can be seen that 104 partial wave contributions are enough to converge the ICS in the collision-energy range of 0.0 eV–1.0 eV.
Through summing the partial waves and the reaction probabilities of different angular quantum numbers, we can obtain the ICS of the titled reaction, which is shown in Fig. 6 for the collision energy increasing up to 1.0 eV. From this figure we can find that the ICS not only shows no threshold due to the exothermal feature of the titled reaction, but also decreases with the collision energy increasing. In addition, the oscillatory character can be found in the low energy region, and it gradually vanishes as the collision energy increases. Figure 6 also shows the experimental results[4] of the ICS (300 K) for the titled reaction. Two blue dashed curves respectively represent the maximum and minimum values of the experimental data, which correspond to a 20% experimental error. We can obviously observe that the present ICS is in good agreement with the experimental results[4] within error margins, and only for
and
are the quantum results in the lower limit of the experimental data.
4. ConclusionsIn this work, we acquire the quantum dynamic information of the O+ +D2 reaction based on a novel accurate global PES of the ground-state H2O+(
) by the Chebyshev wave packet propagation method. Both the reaction probabilities and the ICS of the titled reaction in a collision energy range of 0.0 eV–1.0 eV are provided. The characteristics of the reaction probability curves are in accordance with those of the exothermic reaction with a potential well. Meanwhile it can be seen that the vibrational excitation of the reactant D2 can promote the reaction, and the rotational excitation of the reactant D2 can inhibit the reaction, due to the exothermic and barrierless reaction. The obtained ICSs are in excellent agreement with the existing experimental data, demonstrating the accuracy of the dynamical calculations.
Acknowledgment
The authors thank Prof. Shiying Lin of Shandong University for providing the computational code.
